      SUBROUTINE ANGLAP(SEGSIN,NSEGSIN,SEGSOUT,MAXSEGSOUT,NSEGSOUT,
     *      LUERR,IERR)
      IMPLICIT REAL*8 (A-H,O-Z)
C
C  THIS ROUTINE FINDS THE SECTORS OF A CIRCLE THAT ARE COMMON TO
C  ALL SECTORS INPUT TO IT.
C
C
C VARIABLE  DIM     TYPE  I/O  DESCRIPTION
C --------  ---     ----  ---  -----------
C
C SEGSIN  2,NSEGSIN  R*8   I   THE SECTORS OF A CIRCLE FOR WHICH
C                              REGIONS COMMON TO ALL ARE WANTED.
C
C                              SEGSIN(1,I)=START ANGLE OF THE I'TH
C                              SECTOR.
C
C                              SEGSIN(2,I)=END ANGLE OF THE I'TH
C                              SECTOR.
C
C                              SECTOR(2,I) < SECTOR(1,I) INDICATES
C                              THAT THE SECTOR CROSSES TWOPI.
C
C                              SECTOR(1,I)=SECTOR(2,I) MEANS THE
C                              SECTOR HAS ZERO WIDTH.
C
C                              THE VALID RANGE IS ZERO TO TWOPI,
C                              INCLUDING ZERO AND TWOPI.  IN RADIANS.
C
C NSEGSIN    1      I*4    I   NUMBER OF SECTORS INPUT.  MAY BE ZERO, IN
C                              WHICH CASE NO COMMON SEGMENTS EXIST AND
C                              NSEGSOUT=0 IS RETURNED.
C
C                              MAX=361, BUT CAN BE CHANGED BY REDEFINING
C                              THE MAXIN PARAMETER BELOW AND RECOMPILING
C                              THIS ROUTINE.
C
C SEGSOUT           R*8    O   SECTORS FOUND COMMON TO ALL INPUT SECTORS
C      2,MAXSEGSOUT
C                              SEGSOUT(1,I) IS START ANGLE OF THE
C                              I'TH COMMON SECTOR. SEGSOUT(2,I) IS THE
C                              END ANGLE.  IN RADIANS.
C
C                              SECTORS ARE OUTPUT IN ASCENDING ORDER.
C                              THAT IS SEGSOUT(1,1)<SEGSOUT(2,1)<
C                              SEGSOUT(1,2)<.....<SEGSOUT(2,NSEGSOUT).
C
C                              ALL SEGSOUT VALUES ARE IN THE ZERO TO
C                              TWOPI RANGE, EXCEPT WHEN A COMMON SEGMENT
C                              CROSSES TWOPI(IN WHICH CASE SEGSOUT(1,1)
C                              IS NEGATIVE.
C
C                              SECTORS MAY HAVE ZERO WIDTH. FOR
C                              EXAMPLE, THIS ROUTINE WILL SAY THAT
C                              THE INPUT SECTORS [10 DEG, 20 DEG] AND
C                              [20 DEG, 15 DEG] HAVE TWO COMMON SECTORS,
C                              ONE BEING [10 DEG, 15 DEG] AND THE OTHER
C                              [20 DEG, 20 DEG].
C
C                              THE CALLING PROGRAM CAN USE THE SAME
C                              ARRAY FOR SEGSIN AND SEGSOUT IF DESIRED.
C                              CALL ANGLAP(SEGS,N1,SEGS,N2,N3,K1,K2) IS
C                              VALID.
C
C MAXSEGSOUT 1      I*4    I   MAX NUMBER OF PAIRS THE SEGSOUT ARRAY
C                              CAN HOLD.  IF NSEGSOUT EXCEEDS THIS,
C                              ERROR CONDITION 3 EXISTS.
C
C NSEGSOUT   1      I*4    O   NUMBER OF SEGMENTS FOUND TO BE IN COMMON
C                              WITH ALL INPUT SEGMENTS.  POSSIBLE
C                              VALUES ARE -1,0,1,2,...,NSEGSIN.  THE
C                              NEGATIVE VALUE IS RETURNED ONLY WHEN
C                              AN ERROR OCCURS(IE, IERR NON-ZERO).
C                              
C LUERR      1      I*4    O   FORTRAN UNIT NUMBER FOR ERROR MESSAGES.
C                              = 0 OR NEGATIVE MEANS NO MESSAGE POSSIBLE
C
C IERR       1      I*4    O   ERROR RETURN FLAG.
C                              = 0 MEANS NO ERROR FOUND.
C                              = 1 MEANS NSEGSIN IS TOO LARGE.
C                              = 2 MEANS ZERO TO TWOPI REQUIREMENT
C                                  ON SEGSIN VALUES HAS BEEN VIOLATED.
C                              = 3 MEANS MAXSEGSOUT IS TOO SMALL
C
C***********************************************************************
C
C  C PETRUZZO. 8/82.
C
C  MODIFIED.....
C        12/82. CJP. BUG CORRECTION ASSOCIATED WITH 0, TWOPI ENDPOINTS.
C         9/86. CJP. BUG CORRECTION AND COMMENT MODS TO CLARIFY THE
C                    DESCRIPTION OF THE LOGIC.  MINOR CODE CLEANUP.
C                    PREVIOUS VERSION SAID, ERRONEOUSLY, THAT SEGSOUT
C                    DIMENSION WAS (2,3).  SOME CODE MODS DONE TO
C                    CORRECT THE LOGIC IMPLICATIONS OF THAT ERROR.
C         3/87. CJP. ADDED MAXSEGSOUT ARGUMENT. LARGE INPUT ARRAYS
C                    CAN BE INPUT.  THE MOD ALLOWS SMALL OUTPUT ARRAYS
C                    TO BE USED AND AN ERROR CHECK TO SEE THAT THEY
C                    ARE NOT TOO SMALL.  INCREASED MAXIN VALUE TO
C                    ACCOMMODATE QUIKVIS PROGRAM NEEDS.
C
C***********************************************************************
C
      REAL*8 SEGSIN(2,1),SEGSOUT(2,1)
C
      PARAMETER MAXIN=361,MAXIN1=MAXIN+1
      INTEGER*2 IORD(2*MAXIN)
      REAL*8 SEGIN(2,MAXIN),SEGOUT(2,MAXIN1)
      LOGICAL*1 XTWOPI(MAXIN),LTEST,OKSOFAR,ENDPT,OKTHIS
      REAL*8 TWOPI/ 6.283185307179586D0 /
      REAL*8 DEGRAD/ 57.29577951308232D0 /
      INTEGER NERR/0/
      REAL*8 TOLER/1.D-10/
C
      IERR = 0
      IBUG = 0
      LUBUG = 19
C
      IF(IBUG.NE.0) WRITE(LUBUG,1001) NSEGSIN,
     *    ((SEGSIN(I,J)*DEGRAD,I=1,2),J=1,NSEGSIN)
 1001 FORMAT(/,' ANGLAP. NUMBER OF SEGMENTS INPUT=',I4,' SEGS='/,
     *    99(3(1X,F7.3,1X,F7.3,4X)/))
C
C
C INITIALIZATION
C
      NSEGOUT = 0
C
C    IF NO SEGS WERE INPUT, NONE ARE OUTPUT.
C
      IF(NSEGSIN.LE.0) GO TO 999
C
C ERROR CHECKS
C
C    NUMBER OF SEGS INPUT.
C
      IF(NSEGSIN.GT.MAXIN)  THEN      ! EXCEEDS SUBROUTINE CAPACITY
        IERR = 1
        IF(LUERR.GT.0) CALL MESAGE(1,NERR,50,0,1,0,LUERR,
     *        'ANGLAP. TOO MANY ANGLE SEGMENTS. NO PROCESSING DONE.')
        GO TO 999
        END IF
C
C    RANGE OF VALUES.  [0,TWOPI] ONLY VALID RANGE.  MODULO ARITHMETIC
C    IS NOT DONE HERE BECAUSE THE CALLING PROGRAM SHOULD RESOLVE THE
C    AMBIGUITY TO ENSURE THE PROPER INTERPRETATION OF WHAT START AND
C    END ANGLES MEAN.
C
      DO I=1,NSEGSIN
      DO J=1,2
        TEMP = SEGSIN(J,I)
        IF( TEMP.LT.0.D0 .OR. TEMP.GT.TWOPI+TOLER) THEN
          IERR = 2
          IF(LUERR.GT.0) CALL MESAGE(1,NERR,50,0,1,0,LUERR,
     *      'ANGLAP. INPUT ANGLES OUTSIDE [0,2PI]. NO PROCESSING DONE.')
          GO TO 999
          END IF
        END DO
        END DO
C
C
C  MORE INITIALIZATION
C
C    TRANSFER SEGSIN ARRAY TO INTERNAL ARRAY SEGIN.  MODIFY CONTENTS IF
C    THEY ARE CLOSE TO ZERO OR TWOPI TO MAKE THEM EXACTLY ZERO AND TWOPI
C
      DO J=1,NSEGSIN
        TEMP1 = SEGSIN(1,J)
        TEMP2 = SEGSIN(2,J)
        IF(TEMP1.LE.TOLER) THEN
          TEMP1 = 0.D0
          IF(TEMP2.LE.TOLER) TEMP2 = 0.D0
          IF(TEMP2.GE.TWOPI-TOLER) TEMP2 = TWOPI
        ELSE  IF(TEMP1.GE.TWOPI-TOLER) THEN
          TEMP1 = 0.D0
          IF(TEMP2.LE.TOLER .OR. TEMP2.GE.TWOPI-TOLER) TEMP2 = 0.D0
          END IF
        IF(TEMP2.LE.TOLER) THEN
          TEMP2 = 0.D0
          IF(TEMP1.GE.TOLER)  TEMP2 = TWOPI
          END IF
        SEGIN(1,J) = TEMP1
        SEGIN(2,J) = TEMP2
        END DO
      NSEGIN = NSEGSIN
C
C    DETERMINE WHETHER THE SEGMENTS CROSS TWOPI.
C
      DO I=1,NSEGIN
        XTWOPI(I) = SEGIN(2,I) .LT. SEGIN(1,I)
        END DO
      IF(IBUG.NE.0) WRITE(LUBUG,1003) (I,XTWOPI(I),I=1,NSEGIN)
 1003 FORMAT(/,' TWOPI CROSSING FLAG:'/,99(10(2X,I3,'-',L2)/))
C
C
C
C******* NOW TEST EACH END POINT TO SEE IF IT IS INCLUDED IN ALL OF THE
C        INPUT SEGMENTS.  IF SO, THEN IT IS THE START OR END POINT OF A
C        COMMOM SEGNEMT.  IF NOT, THEN IT IS IGNORED AS A POSSIBLE
C        START/END POINT.
C
C     FIRST, SEE IF THERE ARE ANY ZERO-WIDTH SEGMENTS. IF SO,
C     THIS IS A SPECIAL CASE AND IS HANDLED SEPARATELY. A ZERO-WIDTH
C     SEGMENT, IF COMMON TO ALL OTHER SEGMENTS, IS THE ONLY ONE 
C     POSSIBLE FOR OUTPUT.
C
      ZANGLE = -999.D0
      I = 0
      DO WHILE (ZANGLE.LT.0.D0 .AND. I.LT.NSEGIN)
        I = I+1
        IF( SEGIN(1,I).EQ.SEGIN(2,I) ) ZANGLE = SEGIN(1,I)
        END DO
      IF(ZANGLE.GE.0.D0) THEN    ! A ZERO-WIDTH ONE WAS INPUT.
        DO I=1,NSEGIN
          SEG1 = DMIN1(SEGIN(1,I),SEGIN(2,I))
          SEG2 = DMAX1(SEGIN(1,I),SEGIN(2,I))
          ENDPT = ZANGLE.EQ.SEG1 .OR. ZANGLE.EQ.SEG2   ! AN ENDPOINT?
          LTEST = SEG1.LT.ZANGLE .AND. ZANGLE.LT.SEG2  ! AN INTERIOR PT?
          IF(XTWOPI(I)) LTEST = .NOT.LTEST
          LTEST = LTEST .OR. ENDPT   ! INTERIOR OR ENDPOINT? MUST BE YES.
          IF(.NOT.LTEST) GO TO 999    ! NOT IN ONE OF THE SEGMENTS.
          END DO
C       TO GET HERE, ANGLE VALUE HAD TO BE IN ALL SEGMENTS.
        NSEGOUT = 1
        SEGOUT(1,1) = ZANGLE
        SEGOUT(2,1) = ZANGLE
        GO TO 999
        END IF
C
C  DETERMINE THE ASCENDING ORDER OF THE SET OF ENDPOINTS.
C
      CALL ORDER8(SEGIN,2*NSEGIN,IORD)
      IF(IBUG.NE.0) WRITE(LUBUG,1002)  
     *      (SEGIN(IORD(I),1)*DEGRAD,I=1,2*NSEGIN)
 1002 FORMAT(/,' ENDPOINTS IN ASCENDING ORDER:'/,
     *    1X,F7.3,9(1X,F7.3)/)
C
C  NOW CHECK EACH OF THE END POINTS TO SEE WHETHER IT IS INCLUDED IN
C  EVERY INPUT SEGMENT.  IF SO, SET NSEGOUT AND SEGOUT APPROPRIATELY.
C
      NFOUND = 0
      NSEGOUT = 0
      OLDANG = -1.D10
      NANGS = 2*NSEGIN+2   ! WANT TO TEST 0 AND TWOPI ALSO.
      DO 100 IANG=1,NANGS
C
C     FIRST, DETERMINE THE CURRENT END POINT TO BE TESTED AND SET ISUB
C     TO INDICATE WHETHER IT IS THE START(ISUB=1) OR END(ISUB(2) POINT
C     OF AN INPUT SEGMENT.
      IF(IANG.EQ.1 .OR. IANG.EQ.NANGS) THEN
        IF(IANG.EQ.1) THEN
          ANGLE = 0.D0
          ISUB = 1
        ELSE
          ANGLE = TWOPI
          ISUB = 2
          END IF
      ELSE
        INDX = IORD(IANG-1)
        ISUB = 1
        IF(MOD(INDX,2).EQ.0) ISUB = 2
        ISEG = (INDX+1)/2
        ANGLE = SEGIN(ISUB,ISEG)
        END IF
      IF(IBUG.NE.0) WRITE(LUBUG,1004) IANG,ANGLE*DEGRAD
 1004 FORMAT(/,' TESTING WITH IANG=',I3,' ANGLE=',F7.3)
C
C     CHECK EACH INPUT SEGMENT TO SEE IF THE CURRENT ENDPOINT IS IN ALL
C     OF THEM.  IF SO, IT IS THE BEGINNING OR END OF A COMMON SEGMENT.
      JSEG = 0
      OKSOFAR = .TRUE.
      DO WHILE ( JSEG.LT.NSEGIN .AND. OKSOFAR )
        JSEG = JSEG + 1
        SEG1 = DMIN1(SEGIN(1,JSEG),SEGIN(2,JSEG))
        SEG2 = DMAX1(SEGIN(1,JSEG),SEGIN(2,JSEG))
        ENDPT = ANGLE.EQ.SEG1 .OR. ANGLE.EQ.SEG2  !ENDPOINT ALWAYS 'IN'
        LTEST = SEG1.LT.ANGLE .AND. ANGLE.LT.SEG2
        IF(XTWOPI(JSEG)) LTEST = .NOT.LTEST
        OKSOFAR = LTEST .OR. ENDPT
        IF(IBUG.NE.0) WRITE(LUBUG,1005) ANGLE*DEGRAD,JSEG,SEG1*DEGRAD,
     *      SEG2*DEGRAD,XTWOPI(JSEG),LTEST,ENDPT,OKSOFAR
 1005   FORMAT(/,' ANGLE=',F7.3,' TESTING AGAINST SEG # ',I3,
     *      ' SEG1,SEG2=',2F8.3/,
     *      ' XTWOPI=',L2,'  LTEST=',L2,'  ENDPT=',L2,' OKSOFAR=',L2)
        END DO
C
      OKTHIS = OKSOFAR
C
C     IF OKTHIS=.T., THEN THIS START/END OF AN INPUT SEGMENT IS AN END
C     OR INTERIOR POINT OF ALL OTHER INPUT SEGMENTS.  THEREFORE, IT IS
C     THE START OR END OF A COMMON SEGMENT.
C
      IF(OKTHIS) THEN
C      NFOUND IS THE NUMBER OF UNIQUE ANGLE VALUES COMMON TO ALL INPUT
C      SEGMENTS.
        IF(ANGLE.NE.OLDANG) NFOUND = NFOUND+1
        NSEGOUT = (NFOUND+1)/2
        SEGOUT(ISUB,NSEGOUT) = ANGLE
        IF(IBUG.NE.0) WRITE(LUBUG,1006) ISUB,NSEGOUT,NFOUND
 1006   FORMAT(/,' FOUND ONE. ISUB=',I2,'  NSEGOUT=',I2,' NFOUND=',I2)
        END IF
C
      OLDANG = ANGLE
      IF(DABS(OLDANG-TWOPI).LE.TOLER) OLDANG = TWOPI
  100 CONTINUE
C
C ALL COMMON SEGMENTS HAVE BEEN IDENTIFIED.  IT IS POSSIBLE THAT A
C COMMON SEGMENT CROSSES TWOPI.  IF SO, THEN THE SEGOUT ARRAY REPRESENTS
C IT AS TWO SEGMENTS.  FOR EXAMPLE(IN DEGREES), IF THE COMMON SEGMENT IS
C 300 TO 30(THRU 360), THEN SEGOUT(*,1) = 0,30 AND SEGOUT(*,NSEGOUT) =
C 300,360.  WE WANT TO CHANGE THIS REPRESENTATION TO A SINGLE SEGMENT
C BY REPLACING THE TWO WITH SEGOUT(*,1) = -60,30.
C
      IF(NSEGOUT.GT.1) THEN
        IF( SEGOUT(1,1).EQ.0.D0 .AND. SEGOUT(2,NSEGOUT).EQ.TWOPI) THEN
          SEGOUT(1,1) = SEGOUT(1,NSEGOUT) - TWOPI
          NSEGOUT = NSEGOUT-1
          END IF
        END IF
C
      IF(NSEGOUT.GT.MAXSEGSOUT) THEN  ! NOT ENEOUGH SPACE IN SEGSOUT
        IERR = 3
        IF(LUERR.GT.0) CALL MESAGE(1,NERR,50,0,1,0,LUERR,
     *      'ANGLAP. MORE ANGLE SEGMENTS NEEDED THAN ARRAY HAS SLOTS.')
        GO TO 999
        END IF

C
C
  999 CONTINUE
C
C WRAP UP.
C
      IF(IERR.EQ.0) THEN
        NSEGSOUT = NSEGOUT
        IF(NSEGSOUT.GT.0) THEN
          DO I=1,NSEGSOUT
            SEGSOUT(1,I) = SEGOUT(1,I)
            SEGSOUT(2,I) = SEGOUT(2,I)
            END DO
          END IF
      ELSE
        NSEGSOUT = -1
        END IF
C
C
      IF(IBUG.NE.0) WRITE(LUBUG,1007) IERR,NSEGSOUT,
     *        ((SEGSOUT(I,J)*DEGRAD,I=1,2),J=1,NSEGOUT)
 1007 FORMAT(/,' FINISHED ANGLAP. IERR=',I2,'  NSEGSOUT=',I2/,
     *    ' SEGSOUT=',(T12,2X,2F8.3,2X,2F8.3,2X,2F8.3))
      RETURN
      END
